Analytic static structure factors and pair— correlation functions 
for the unpolarized homogeneous electron gas 
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We propose a simple and accurate model for the electron static structure factors (and correspond- 
ing pair-correlation functions) of the 3D unpolarized homogeneous electron gas. Our spin-resolved 
pair-correlation function is built up with a combination of analytic constraints and fitting proce- 
dures to quantum Monte Carlo data, and, in comparison to previous attempts (i) fulfills more known 
integral and differential properties of the exact pair-correlation function, (ii) is analytic both in real 
and in reciprocal space, and (iii) accurately interpolates the newest, extensive diffusion-Monte Carlo 
data of Ortiz, Harris and Ballone [Phys. Rev. Lett. 82, 5317 (1999)]. This can be of interest for 
the study of electron correlations of real materials and for the construction of new exchange and 
correlation energy density functionals. 



I. INTRODUCTION 

The homogeneous electron gas, a model solid whose 
positive ionic charges are smeared throughout the whole 
crystal volume to yield a shapeless, uniform positive 
background (whence the nickname of jellium) has pro- 
vided, since the very start of quantum mechanics, a key 
conceptual reference and a mine, of information for solid- 
state and many-body theorists. ETu Initially it was mostly 
regarded as an approximation of the true distribution of 
valence electrons in simple metals, since, in spite of its 
crudity, it could already account for some of their experi- 
mental properties .□ Although the importance of valence- 
charge inhomogeneites in real materials wasL soon rec- 
ognized (and described first by perturbationEl and later 
by self-consistent pseudopotential theoryQ) , the homoge- 
neous electron gas stood by itself, ovex the decades, as. an. 
independent active field of theoretical] and numericalErE£l 
investigation. One reason for this continued interest is 
that the model, by ignoring the ionic lattice which makes 
real materials different from one another, allows the the- 
orists to concentrate on key aspects of the electron- 
electron interaction. Another reason for caring about 
such an unrealistic system residesdn its connection to the 
inhomogeneous electron gas:Ei&tJ not only does the jel- 
lium model represent an obvious-limit, but also, through 
the Density Functional TheoryEj and its Local Density 
Approximation, it links to ajaopular and very successful 
description of real materials S3 For the rattes-reason, from 
the simplest Hartree-Fock approximation^ to the pio- 
neering QMC simulations Jj almost any theory of jellium, 
its electron correlations and its pair-correlation functions 
has also implied an improved understandipg.-and con- 
struction of Kohn-Sham energy functionals .lHHi3 

In this context our work aims at a new simple analytic 
expression for the pair-correlation function of the homo- 
geneous electron gas, which describes the spatial correla- 
tions of electron pairs with prescribed spin orientations. 



A good model pair-correlation function and static struc- 
ture factor has its own interest; its availability over a wide 
density range is crucial for new developments and appli- 
cations of the Density Functional Theory, through the 
construction of ab initio exchange and correlation energy 
functionals in generalized gradi ent approximationsEj and 
in other beyond-LDA schemes. E3E3~E3 As a consequence, 
over the last 20 years, several authors already pr 
ingenious expressions for this or related functionsEIrET 
A first motivation for resuming and improving over previ- 
ous efforts is the avalaibili ty ., fr .om recent quantum Monte 
Carlo (QMC) simulations ,I12h13 of a wealth of new numer- 
ical results for the pair-correlation functions and static 
structure factors of jellium. A second motivation comes 
from the observation that most of the previous models 
were not spin-resolved and all of them neither fulfilled all 
the known exact properties, nor were given in analytic, 
closed form both in real and reciprocal space. Our goal is 
thus to give a new, spin-resolved expression for the pair- 
correlation function which is analytic both in real and 
reciprocal space, automatically incorporates more exact 
properties than any previous expression, and contains 
enough free parameters to fit the new QMC resultsO 

We recall the exact properties of the pair-correlation 
function for the unpolarized jellium in Sec. [ll| The three 
subsequent sections are devoted to a description of our 
general strategy (Sec. Ill) and of the resulting functional 
form for the antiparallel- (Sec. [[^) and parallel-spin 

VJ we de- 
Oncc the 



(Sec. ^) pair-correlation functions. In Sc< 
scribe our fitting procedure to QMC data 
exact constraints are imposed, 18 free parameters (9 for 
antiparallel spins and 9 for parallel spins) are enough 
to yield extremely accurate two-dimensional fits of the 
w 9000 + 9000 new QMC data pointsB as a function of 
the interelectronic distance r and the density parameter 
r s in the relevant density range r s < 10. Our results are 
discussed and compar ed w ith the widely used Perdew- 
WanepJ model in Sec. VII. 
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In Sec. VII] we report the correlation energy which 
corresponds to our model pair-correlation function-n and 
we find that its agreement with the QMC energies)^ not 
targeted by our fitting procedure, is as good (~ 5%) as 
the most popular interpolation formulae for the correla- 
tion energy. The last Sec. IX is devoted to conclusions 
and perspectives. 



g(p\r s ) = - [g^(p;r s ) + gn(p;r s )] 
S(k; r s ) = S n (k; r s ) + S n (k; r s ). 



B. Pair—correlation function near r = 



(5) 
(0) 



II. EXACT PROPERTIES 

We briefly recall many known properties of the pair 
correlation function of the unpolarized homogeneous elec- 
tron gas. Its integral properties (sum rules) will be 
rewritten in terms of q — > properties of Fourier trans- 
forms, since this choice turns out to be convenient for our 
subsequent steps. Hartree atomic units are used through- 
out this work. 



A. Definitions 

For an electronic system the pair-correlation function 
g (T1(T2 (r 1 , r 2 ), if n a (r) is the density of electrons with spin 
a =t or I, is defined by 

n ai (r 1 )n r72 (r2)ga 1 a 2 (ri,r 2 ) = 

($1^(^)^(^)^3(1-2)^ (ri)|*> (1) 

and is thus related to the probability of finding two elec- 
trons of prescribed spin orientations at positions i"i and 
r 2 . The normalization of g is such that the expected 
number of electrons of spin a 2 in the volume dV at r 2 
when another electron of spin <ii is known to be at ri is 
given by 



dN(v 2 <7 2 \rii7i) 



(r 2 )g <71< r 2 (r 1 ,r 2 )dV; 



(2) 



the lack of any correlation amounts, then, to the con- 
dition g aia - 2 (ri,r 2 ) = 1. In the spin-unpolarized jel- 
lium the electronic spin density rif(r) = ri|(r) = n/2~ 
(87rrf/3) _1 is uniform in space (i.e. independent of r),E3 
so g ai a 2 ( r i ; r 2) only depends on the distance between the 
two electrons r = |ri — r 2 \. The static structure factor- 
s'^) is directly related to the Fourier transform of the 
pair-correlation function. For an unpolarized homoge- 
neous electron gas, after introducing the Fermi wavevec- 
tor q F = (3TT 2 n)^ 3 = a/r 8 , with a = (97T/4) 1 / 3 , the 
scaled variables p = q F r and k = q/q F are often conve- 
nient. With these variables the static structure factors 
are written as 



S n (k) 



2 

2 

37 



1 



d P \9m{p) - 1] P 
dp [9TT(P) - !] P 



2 sin (M 



k p 

2 s m(kp) 
k p 



(3) 
(4) 



and the total pair-correlation function and static struc- 
ture factor are given by: 



The behavior of g aia2 (r;r s ) in the r — > limit can 
be directly obtained from the many-body Schrodinger 
equation whea two electrons approach each other (cusp 
conditions) :E3H 



g;9n\r,r B ) 
d 3 



= flu(r->0;r.) ^0 



= <?Tt( r °; r s) = 



3 ^ ( \ 



(7) 
(8) 
(9) 



Eqs. (0)-(||) hold for any 3D-system of N fermions in- 
teracting via the two-body repulsive Coulomb potential. 



C. Structure factor near q — 

The conservation of particles in the system implies the 
relations: 



Su(?->0;r.)=£ n ( ? -0;r 4 ) = 0. 



(10) 



The asymmetry between the definitions (|^) and (||) leads 
to the two well-known sum rules for g^ and (see for 
instance Ref. p3| ). 

The long-wavelength behavior of the total static struc- 
ture factor of Eq. (|^) is determined by the plasmon 
contribution, proportional to q 2 , and by the single-pair 
and multipair quasiparticle-quasihole excitatiqxL_contri- 
butions, proportional to q 5 and q 4 respectivelyffij 



S(q^0;r s 



2u> p (r a ) 



+ Cq 4 + 0(q 5 ), 



(11) 



where to p (r s ) = y3/Vf is the classical plasma frequency. 
In the paramagnetic gas, the parallel and antiparallel- 
spin contribution to the plasma mode is the same. More- 
over, to build up model functions for the spin-resolved 
S ai<72 , it's crucial to include the following property of the 
so-called magnetic structure factor 5ff — S^y. 



S 



n 



-k 



16 



0(k% 



(12) 



where the scaled variable k — q/q F has been used. 
Eq. (113) is valid in the framework of the random-phase 
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approximation (RPA) and can be obtained from a se- 
ries expansion of Sf\ (k) near fc = (Ref. [56|), and 



u 

from the corresponding expansion of the total 5' RPA (see 
for instance rfief. |3^) . Since in the k — > limit RPA 
is exact ,Hl23~Ea we expect Eq. (|l2|) to hold for the exact 
structure factor as well. From Eqs. ([fl]) and ([l2|) we can 
write the small-fc expansion of S aia2 : 



S 



11 



k — 



= lk + 



4w p (r s ) 

q 2 F > 

4wp(r s ) 



32 



0(fc 4 ) 



32 



(13) 



(14) 



D. Correlation energy 

The electron-electron potential energy is, as known, 
given by the sum of repulsive two-body Coulomb terms: 



N 



u = 



-y- 



i 



(15) 



Its ground-state expectation value (per electron), in a 
homogenous electron gas of density r s , is given by the 
following integral over the pair-correlation function: 



(U)r s = 



[g(p;r s ) - 1] pdp. 



(16) 



2 a 2 r s Jo 

By the virial theorenS and the usual definition of the 
correlation energy e c as total electronic energy minus 
Hartree-Fock energy, we have: 



Mr.=-^ + ~[r>ec(r.)] 
4 7r r, ar. 



(17) 



Putting together Eqs. (|16|) and (|17j) one obtains the 
known relation between g(p; r s ) and the exchange and 
correlation energy: 



3gF 
47T 



2a 2 r? 



dpp I dr' s [g(p\ r' s ] 
la 



-1] 
(18) 



The same relation can be obtained in a more gen- 
eral wayt-3 by the Hellmann-Feynman theorem and the 
coupling-constant average of g(p; r s ), which, for the ho- 
mogeneous system is just the average over r s : 



g(p; r s ) = — / g(p; r' s )dr' s 



indirectly related to the exchange 
Ej of the electron gas. 



The function ~g(p\ r Jl 
and correlation holdllr 

We have recalled these relations because we will later 
check our analytical expressions for g(p; r s ) against avail- 
able energy data, and also because, among other con- 
straints, we want our functional form of g(p;r s ) to be 
consistent with the high-density limit of e c (r s ): 



e c (r s -vO) -A In r s + B + Cr s \nr s + Dr s (20) 

where A, B, C and D are known constants, 

A = (1 - ln2)/vr 2 , B = -0.0469205, C = 0.0092292, 

D = —0.01, and the next leading term is 0(r 2 s lnr s ). 



III. GENERAL STRATEGY 



We study the antiparallel and parallel-spin correlation 
functions both in real and reciprocal space and we split 
them, as usually, into exchange and correlation according 
to: 



9ll(p;r„) = l+g c tl (p;r s ) 

9n(P: r s) = 9ex(p)+g^(p;r s ) 
s n( k ' r s) = Sn(k;r s ) 
S n (k;r s ) = S ex (k) + Sfr(k;r„) 



(21) 
(22) 
(23) 
(24) 



where the exchange functions, given by the Hartree-Fock 
approximation, are equal to: 



>(p)=l-9 



sin p — p cos p 



5 e ,(fc) = | 3fc / 4 - fc3 / 16 forfc ^ 2 



[ 1 



for k > 2 



(25) 
(26) 



and our model only concerns the correlation part. 
Putting together Eqs. (§), @, @, @ and @ one 
finds a well-known result: in the total S = S ex + + 
the linear term of S ex (k), 3k/ 4, which dominates its 
small-fc behavior (and corresponds to a large-p leading 
term cx 1/p 4 of g ex ) exactly cancels the small-fc lead- 
ing term of the correlation part SS^ + SZ+. This prop- 
erty has been incorporated in several previous functional 
forms for the total g, as the widely used Perdew-Wango 
model (hereafter PW) where, however, the k 2 coefficient 
in the small-fc expansion of S(k) is slightly different from 
the exact one [Eq. (O)] , because of spurious fc 2 contribu- 
tions from their (g^)[Eq. (19) of Ref. p9fl and from their 
short-range part of g c [Eq. (37) of Ref. [29|. 

The fc -> limit of Eq. @ seems, instead, to-.be less 
known: even the best-to-date spin-resolved PWE^I model 
does not incorporate such a non-trivial analytic prop- 
(19) erty, which can alternatively be expressed as SS^ beeing 
identical to in the the small-fc limit and corresponds 
to a visible feature of the magnetic structure factor (see 
Sec. VII). Our goal is to produce simple and practical 



analytical functional forms for SS,(k;r s ) and S°+(k;r s ) 
[and hence g^Ap, r s ) and g^j{p,r s )] which satisfy all the 
physical properties of Sec. O and have enough variational 
flexibility to accurately interpolate the QMC data of Or- 
tiz, Harris and Ballonc.ll3 



3 



To do this, let's start with a few simple considera- 
tions about spherical Fourier transforms: S° 1<T2 (k) and 
9a 1 a 2 jj ) ) are related to one another by an integration like 
Eq. (j3J). The function sin(fcp)/fcp is an even function, i. 
e. its odd derivatives in fc = (or p — 0) are all equal to 
zero. However, the small-/? and the small-fc properties of 
g and S tell us that they must have non-zero odd deriva- 
tives in p = and k — 0. This is achieved if (and only 
if), as the integration variable goes to infinity, the inte- 
grand goes to zero slowly enough as to avoid absolute 
convergence, so that differentiation within the integral 
sign is not allowed. It is easy to establish a connection 
between the large-fc (large-/?) behavior of S (g) and the 
odd derivatives in p = (k = 0) of g (S): a derivative of 
g in p = of order 2n + 1 corresponds, in S, to a large- 
k term oc l/fc 2ra+4 and viceversa. This simple relation 
was used in Ref. ^3] to obtain the large-fc expansion of 
S aia2 from the cusp conditions of Eqs. (Q) and ([)]). These 
elementary considerations lead us to write down a very 
simple functional form S c in reciprocal space which au- 
tomatically has the exact small-fc and large-fc behavior. 
Its spherical Fourier transform g c is analytic and closed- 
form, consisting of the same kind of functions used in 
reciprocal space. We thus have an equally simple expres- 
sion for S c and g c . 

We begin by studying the antiparallel-spin part, and 
do it in several ste ps. F irst (IVA) we choose our func- 



n 



tional form. Then ( IV B ) we impose to it the properties 
of Sec. y. At this point we are left with 6 free parame- 
ters, which, independently for each available r s , are used 
to accurately fit the QMC data both in real and recipro- 
cal space, as done in Ref. [h]. In our case, however, the 
r s dependence of each of the 6 optimal parameters turns 
out to be both regular and monotonic. We then try to 
represent each of them as a simple function of r s in such a 
way that (i) as r s — > the exact high-density exp ansio n 
of the correlation energy [Eq. (20)] is recovered ( [ V C ) , 
and (ii) for—finite r s < 10 an optimal global fit of all the 
QMC dataO is obtained ([vf). We apply the same strat- 
egy to the parallel-spin part (|y|). Besides the excellent 
quality of the final fits of g and S, we see that even the 
resulting correlation energy, not targeted by our fitting 
strategy except at r s — ► 0, turns out to be in good agree^ 
ment (within 5%) with the corresponding QMC resultsEa 
at any r s . We compare our correlation energy with the 
most popular interpolation formulae and we discusSr-their 
relative efficiency in fitting the new QMC energies.ll3 



IV. ANTIPARALLEL SPINS 

A. Functional form 

In reciprocal space our functional form is simply writ- 
ten as: 

6 

S^(k; r s ) = exp [~b n (r s ) k] ^ k n + 

n=l 



(r.) k 8 



(r s )k 



in 



[(aU)2 + fc2] 



(27) 



as mentioned, the corresponding gi^ amounts to a linear 
combination of the same kind of functions0 in real space 



[see Apr 
in Eq. 027H 



Eq. (Al)]. Two types of functions appear 
the first one, an exponential cut-off times 
a truncated power series, fully characterizes the long- 
wavelength behavior of S, while the second one entirely 
determines its large-fc expansion. The leading term as 
k — * oo isjaf-order fc -4 , as exactly known from the cusp 
conditionjEj'El in real space the short-range behavior is 
thus entirely determined by the parameter a\^(r s ): 



d_ 

dp 



9u(p = °; r s) 



3tt 

T 

3?r 
T 



a\ l (r s ) (28) 
q F (r s )al L (r s ). (29) 



B. Physical constraints 



The fc — > conditions of Subsec. [I C are easily im- 
posed: 



= b n c n + Qf_ 



4w<o 



4w p 32 



(30) 
(31) 
(32) 



The cusp condition of Eq. (m) fixes a simple relation be- 
tween ctg and the other parameters: 



„TJ- - ^il 




(33) 



After imposing all the |J. physical conditions, our model 
[Eq. (27)] is left with 6 free parameters: the two expo- 
nential cut-offs (a^ in real space and in reciprocal 
space), the parameter , which determines the short- 

n 



range behavior of g(r), and the 3 linear parameters c 4 



^ 5 and , which will be used to further increase the 
variational flexibility and, fit the numerical gn obtained 
from QMC simulations.!^ The dependence of these free 
parameters on r s will be determined according to the 
strategy summ arized in Sec. p| and detailed in the fol- 
lowing Subsec. IV C and Sec. lW 
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C. High—density expansion 



As anticipated in Subsec. II D , we want our pair- 
correlation function such that its f J, and ft contributions 
automatically fulfill the high-density limit of the corre- 
lation energy. We thus fix the r s — > limit of our free 
parameters by means of Eqs. Jl^), ( |l7| ) and fl2C|). Our 
antiparallel contribution to the expectation value of the 
potential energy U = Un +U-[-[ is simply given by: 



is worthwhile to point out that this contradiction is due 
to the exponential cut-off times a truncated power series 
in real space [Eq. (|A1|)], and emerges when the cusp con- 
dition of Eq. (Q) is imposed to it. The relevance of this 
limitation, which only c oncer ns densities r s < 0.1, will 
be discussed in next Sec. 



VII. 



V. PARALLEL SPINS 



u n = — 



6 

' /F E 



+ 



2048 



7a 



21a 



11 



(flTi)5 ( a U)3 



(34) 



In the high-density limit, the correlation-energy con- 
straint of Eq. ( p0| ) translates into the following condition 

{Un)r.^0=2A n \nr s + (A n +2B n ) + 0(r s lnr s ) 

(35) 

where, comparing to Eq. (|2^), is simply givenEi by 
\A. To determine B^ we recall that the constant B 
in Eq. ( |2C| ) is the sum of two contributions: a second- 
order exchange term, Be-cc, which only concerns the H 
part, and a direct term, Bd, which is, instead, equally 
split (in the unpolarized gas) between H and ||- Hence, 



exactlyc 



(2) 

Both Bex C and B^ have been evaluated 



Provided that the 3 linear parameters cj , Cr and c, 



remain finite as r s — ► 0, the exact high-density limit of 
Eq. (p5|) amounts to the following conditions: 



TJ, n ^ _ 1 + fc i^ lnr s + k 2 r s + Q(r 2 s lnr s ) 
" 4 ^ Uj " -3*^/4 

^(r s -o) = (A) V3 7r y| +0(r o ) 

a n (r s -> 0) = const. + 0(r s ) = a n + 0(r s ) 



(36) 

(37) 
(38) 



A. Functional form 

For the correlation part of the ft pair-distribution 
function we apply the same strategy used for the 
antiparallel-spin case. In reciprocal space we thus have: 

6 

S^(k;r s ) = exp [-6 TT (r s ) k] £ (r s )k n + 



a\l(r s )k s + aj\r s )k w + a^ T (r s )fc 12 

[( a TT)2 + fc 2]9 



(41) 



which again corresponds, in real space, te a linear com- 
bination of the same kind of functionscZl [see App. [A|, 
Eq. ( |A2|) ]. The long-wavelength term has the same form 
as the tl part. The large-A: term describes the short- 
range behavior of g^y. the cusp condition of Eq. (||) tells 
us that, as k — > oo, the leading term of must be of or- 
der k~ 6 . With respect to the || case, one more parameter 
is needed for the large-fc term to satisfy the Pauli prin- 
ciple. As in the antiparallel-spin case, the short-range 
properties of g^ are characterized by the parameter: 



Q3 

dp 

dp 



2 stt(p; r s) 



3?r 



p^O 



p^O 



q F (r s )al 1 (r s ). 



(42) 
(43) 



where k\ and ki depend on Afj, B^ and a^: 
18vr(aTi) 2 A 

h = — - — -A n 



ko = 



a 

729 (qt+) 2 
64 a 4 



21 1 

64 a^a 



9(a 



m2 7 



2 a 



-(A n +2B n ). 



(39) 



(40) 



B. Physical constraints 

The small-A properties imply, for the |T case, iden- 
tical constraints as for the f J. case [see Eqs. (|30|), ( |3l| ) 
and @]. 

The Pauli principle and the cusp condition of Eqs. (||) 
and (^|) fix the dependence of and a}J on the remain- 
ing parameters: 



Inserting Eq. ( p6| ) into Eq. (29) we see that in our model, 
once the high-density expansion of e c (r s ) is fixed, the 
r s — > limit of g(p = 0;r s ) is also fixed to the larm 
1 + Cr s lnr s + 0(r s ). The corresponding exact forrrO is, 
up to orders r 2 , 1 + C\r s + C^r 2 lnr s , thus slightly dif- 
ferent from ours. Evidently, the simple functional form 
of Eq. ( p7j ) does not correctly describe the short-range 
behavior of the Coulomb hole at very high densities. It 



„TT 



J I 



2048 
3tt 



(« TT ) 5 E 



^ (& TT ) 

n— 1 v ' 



n+3 



(n + 2)! 



5(n + 4)! 
(aTTfoTT)2 



4096 TT3 _ „ n 
33tt v ' 6 y ' \ 33 



l|. 3 / 2560<7f 26 



,,11 



2048 
3tt 



(« TT ) 7 E 



! (& TT ) 



n+3 



(44) 



(«TT & TT)2 3 ^ + ^- 



5 



4096 
' 15vr 



(a TT ) 5 + 



IT 
3 



(a TT)5 ^ + 512(?F 



(45) 



The cusp condition of Eq. (Q) is not included in the PW 
model. As in the antiparallel-spin case, we have 6 free 
parameters: the exponential cut-off in real space, 



,TT ) 
TT 



separately for the Tl and the ft parts. The r s depen- 
dence of the parameters turns out to be quite smooth 
and monotonic and well described by the following func- 
tional forms ( which also take into account the exact 
high-density limit of Eqs. and @ and guar- 

antee the exact lo w-density expansion of the resulting 



the exponential cut-off in reciprocal space, , the a 6 correlation energy 
parameter, which determines the short-range behavior 
of gf-f-f, and the three linear parameters c\\ clj and 
which arc used to fit the oscillatory behavior of g^. 
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-3/2, 



„tt 



C. High— density expansion 

The contribution to the expectation value of the po- 
tential energy due to the correlation part of our 
(?7 TT = -3g F /47r + C/ T c T ) is: 



TT ^ (&TT)»+ 

71— 1 V J 



<1F 

65536 



35 a 



n 



45 a 



11 



(a 



11)9 



(aTT)7 



99 a 



'1 



(aTt)5 



(46) 



4 l-fci(a^)r s ln(l + fc 2 (a^)/r s ) 



37rq F (1 + k 3 r 2 s ) 
[l- Pl (a^)r s ln(l+p2(o^)/r.)] 



5nq F (l+p 3 rf) 



l (r s ) = a CTl£ 



b a ^(r s ) 



4 

9tt 



1/3 



1 



3/2 



n = 4, 5, 6 



(48) 

(49) 
(50) 

(51) 
(52) 



where fci(a^) is given by Eq. (39), and ^2(0^), pi(a^) 
and ^2(0^) are equal to: 



As r s — ► 0, the condition on US^ is identical to Eq. ( |35[ ) 
whereE] A 



hA, and B TT = B { eX 



hBd. As in the 



TT - 2 J 

t| case, the exact r s — > expansion of US* implies for 
the two exponential cut-offs in real (a^) and reciprocal 
space (b^) identical conditions as Eqs. ( ]37| ) and (|38|). For 



Og the condition is similar to Eq. (]36|): 



0) 



4 



- +pir s \nr s + p 2 r s + 



-0{r s \nr s 



(47) 



&2(g ) = exp 



Pi(« TT ) 



B 1i _ ]_ 
An 2J 

337rA TT (a") 4 



P2(a TT ) = exp 



7 



81 



3847r(aTi)3^ u l287ra 3 A 



11 



81 



9607r(aTT)M n 1287ra 3 A TT 



B 11 1 



.4 



11 



(53) 
(54) 

(55) 



where p\ and pi depend on A^, Bjt and a'' through 
equations similar to Eqs. ( |39| ) and (p0|) . From Eqs. ( [43] ) 
and ( |47| ) one can see that, as expected, when r s — ► 0, 
.9||(p — 0) g° es to the Hartree-Fock value 2/5. As 
in the antiparallel-spin case, the high-density limit of 
the correlation energy fixes the r s — > expansion of 
9r\{p — 0;r s ). The exact g"^(p — 0;r s — > 0) should 
have the formEj 2/5 + 0(r s ), while our functional form 
gives 2/5 + 0(r s lnr s ). Again, we find that in real space 
the simple e xpo nential cut-off times a truncated power 
series [Eq. (A2)] does not correctly describe the short 
range Coulomb interactions at very high densities. 



VI. FIT TO QMC DATA 

For each available density in the range 0.8 < r s < 10 
(i.e. r s = 0.8, 1, 2, 3, 4, 5, 8 and 10) we performed, 
a best fit of the 6 free parameters to the QMC data,E2l 



The 9 constants for || and the 9 constants for j"f have 
been fixed by a two-dimensional best fit to the QMC data 
in real and reciprocal space (9368+9368 data points). 
The efficiency of our interpolation scheme has been tested 
by performing preliminary fits in wich some of the avail- 
able r s were not included and then verifying that the 
interpolated g and S were in good agreement with the 
corresponding QMC quantities for the excluded r s . Since 
this was always the case, we included all the available r s 
in order to have optimal values for our final parameters. 
We thus expect our g and S to be very reliable and ac- 
curate in the whole density range r s G [0.8,10]. The 
optimal 9 parameters which define our best antiparallel- 
spin model are: a n = 0.838, k 3 = 0.141, b\ l = 3.27, 



\\ l = -78, 
itl 



u 

74 



= 28, 



\U _ 

A 5 - 



216, 7 5 U = "124, 
Xg l = —140, 7^ = 55, and the other 9 parameters for the 



parallel-spin part are: o" 
\TT _ na „,TT _ q« \TT 

TT 
7e 



1.32, p 3 = 0.015, 6} T =3.47, 



- 98, 7 4 M 
= -13. 



= -36, A 5 = 



295, 7 5 TT =74, A<U =170 



T1 
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VII. RESULTS IN REAL AND RECIPROCAL 
SPACE 



After fixing the 9+9 parameters which fully specify our 
model, we are now ready to present, in Fig. [j], our real- 
space pair-correlation function g aia - 2 , shown as a solid 
line, as a function of the scaled variable r/r s . r^Chis is done 
for the 8 values of r s for which QMC resultsE^I, shown as 
solid dots, were available. The hest-to-date model corre- 
lation function of Perdew-Wangcj (PW) is also shown for 
comparison as a dashed line. Perdew and Wangcj inter- 
polated the total pair-correlation function g(r) between 
its short-range limit, dominated by the on-top value and 
cusp, and the nonoscillatory part of its long-range limit. 
Their interpolation, controlled by normalization and en- 
ergyJntegrals, agreed with older spin-unresolved QMC 
datala. They only needed the total g(r) foe-construction of 
the generalized gradient approximation,!!^ however, they 
also made an estimate for the spin resolution of g, using 
scaling relations that preserve the normalization integrals 
but are exact only for the exchange contribution. 

Our new expression, explicitly constructed to fit spin- 
resolved numerical correlation functions, follows the 
QMC dataliia better (low r s ) or much better (medium 
and high r s ) than the corresponding PW model, .sdiose 
performance with respect to the new QMC dataEl be- 
comes reasonable only after summing the two contribu- 
tions and going back to the total, spin-unresolved ver- 
sion (not shown) . This can be guessed from the fact that 
for r s > 2, where the discrepancies become clearly visi- 
ble, they generally have opposite signs: both the up-up 
and the up-down correlations are larger (i. e. less close 
to one) than they should. This is due to the fact that 
the PW estimate for the |T P ar t is a simple rescaling 
the pair-correlation function of the fully polarized gas 
while in the unpolarized case correlations are dominated 
by | J, interactions (see for instance Ref. ||). Like the 
PW model, our pair-correlation function breaks down 
for r s > 10: for very low densities g tends to become 
negative at small p. This is probably due to the limited 
variational flexibility of the model, which in this low- 
density regime cannot at the same time fulfill the cusp 
conditions at p = and reproduce a flatter and flatter, 
yet non- negat ive g for p > 0. As we shall see in the 
next Sec. VII] , such a breakdown has no impact on the 
resulting correlation energy, which is an integral of g and 
remains accurate at any r s . 

We compar£-|in Fig. || our gn(p = 0;r s ) (solid line) to 
the YasuharaE3 electron-electron ladder approximation 
(dashed line). Built up to interpolate the QMC data, 
our gn(p — 0;r s ) is larger than the Yasuhara result for 
r s > 0.5, as expected from Fig. |l|, where the discrep- 
ancy between the short-range behavior of the QMC data 
and the PW model (which by construction follows the 
Yasuhara approximation) is clearly visible. In the inset 
the corresponding high-depsity expansions are shown, to- 
gether with the exact limited (dots), which, as anticipated 



in Subsec. IV C| , is not fulfilled by our g^[(p = 0;r s ). 
Rather than giving up our exact limit of e c (r s ) for r s — > 
or trying to fulfill both the e c and the 5tj.(0) limits, we 
have preferred to accept a slight discrepancy of g(0) and 
keep our functional form as described up to now: in our 
experience the collateral complications, at least within 
our functional form, were not worth the effort. Because 
of this limitation, our pair-correlation function does not 
fulfill the high-density limit of g c /r s recently computed 
by Rassolov et al£3 

In Fig. H we report the total static structure factor 
Sff + Sn and the magnetic structure factor Sff — for 
the same 8 values of r s as in Fig. [jj jAgain, our model 
is shown as-solid lines, the QMC dataEj as dots and the 
PW modeled as dashed lines. Our combination of ana- 
lytic constraints and fitting procedure nicely interpolates 
the QMC data, filtering out their noise. In reciprocal 
space it becomes clear that the long-range of the PW 
spin-resolved model is not exact. Moreover, as said in 
Sec. Ill, the PW total static structure factor does not re- 
cover, as q — > 0, the exact plasma frequency in its leading 
q 2 /2oj p term. This is visible for r s — 8 and 10. 



VIII. CORRELATION ENERGY 



A. Spin— unresolved 

The correlation energy obtained by integrating our g 
[see Eq. (|l^)1 is reported-in Fig |J, together with the 
corresponding QMC dataJiia Its ft and f [ contributions 
(e c = e]} +eJ T ) are also separately shown. As expected,El 
correlations are dominated by fl interactions. Our to- 
tal correlation energies are in agreement with QMC data 
within 5% (the maximum absolute error is 3.4 mRy). 
Notice that, even if our model pair-distribution function 
breaks down for r s > 10, it gives very good correlation 
energies even at higher r s values. This is due to the op- 
timal choice of the r s dependence of our free parameters, 
which also includes the low-density expansion of e c . 

To have an idea of the accuracy of our correlation 
energies, we performed beat fits of the QMC data of 
Ortiz, Harris and BalloneEa (hereafter OHB) based on 
other popular interpolation formulae for e c (r s ), i. 
the Perdew-ZungerB (PZ), the Vosko-Wilk-Nusainla 
(VWN) and the Perdew-Wangkll (PW2, to distinguish it 
from the pair-correlation model) functional forms. The 
new QMC data for the correlation energy of the unpo- 
larized jellium are available for a large set of r s : 0.8, 1, 
2, 3, 4, 5, 8, 10, 20, 30, 40, 50 and 60. The results are 
the followings: with the PZ formula one obtains a rather 
good fit (within 3%) , but a wrong negative coefficient for 
the high-density term r s In r s , an unpleasant feature al- 
ready pointed out in Ref. [l(]. Moreover, the PZ energy 
has a discontinuity in its second derivative at r s = 1, an 
unpleasant feature for whoever is interested in the cor- 
responding pair-correlation function, related to the first 
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derivative of e c . 

The VWN form efficiently interpolates the OHB data 
(2.7% maximum relative error; 1.5 mRy maximum ab- 
solute error) only if the free parameter xo of the VWN 
formula has a positive value, which however implies an 
unphysical logarithmic divergence at finite r s (~ 0.6). If 
xo is constrained to be negative, then the fit provided by 
the VWN form is not better than ours (5.2% maximum 
relative error; 3.4 mRy maximum absolute error). 

The fit accomplished with the PW2 form is not very 
accurate (see also Ref. |l0): 7% maximum relative er- 
ror, 3.4 mRy maximum absolute error. Moreover, the 
optimal fit parameter (3 4 of PW2 form turns out to be 
negative (see also Ref. |o|), thus leading to a negative co- 

3/2 

efficient for the low-density expansiottuterm r s and to 
the violation of the Ferrel conditions 

We can conclude that the correlation energies vs. r s , 
which directly emerge from our pair-correlation func- 
tions, although not targeted by our fits, are rather good. 

The main inaccuracies of the popular correlation- 
energy models just reviewed are located in the high- 
density— .pgion, where at first sight the new QMC 
resultal3*E2l cannot be reconciled with the exact r s — > 
limiting behavior. This discrepancy can be related to 
the combined impact of fixed-node approximatior£j and 
infinite-size extrapolation (which would match the find- 
ing that, for r s < 2, Monte Carlo simulations based on 
different nodes and size-scaling rules Refs. || and |l2| ob- 
tain somewhat different energies); it should be kept in 
mind, however, that the exact high-density expansion 
only holds for r s — > 0, and could, in principle, start dom- 
inating the correlation energy at smaller r s values than 
implicitly assumed by the existing models. 

An alternative correlation-energy model, not related 
to our pair-correlation function but capable of an ex- 
cellent interpolation of the QMC energies of Refs. [l(] 
and |l3] including those at high density, can be obtained 
by a minor generalization of the PW2 form. Such a 
generalization keeps its exact r s — * limit, improves 
some of its original analytic properties, and appears flex- 
ible enougl^tXLjjjterpolate different sets of high-density 
QMC data.EMO'Ej We separately present it in our Ap- 
pendix H. 



In Fig. |5| we compare our parallel-spin part of the cor- 
relation energy with two copesponding widely-used scal- 
ing guesses: Perdew-Wangea [eJT(r s ,C = 0) — e c(?-s>Cr] 
l)/2 1 / 3 , where < = |ra T - nl\/n] and Stoll et al£3 
[eJ T (r s ,C = 0) = e c (2 1 / 3 r s ,C = 1)]. Both seem to over- 
estimate the ft contribution to the correlation energy. 
Even if the Stall et alx3 estimate fulfills the exact high- 
density limito (A/2 In 7-3), the PWc£l model (in which 
the r s -t limit is violated) seems to do better in the 
relevant density range r s > 0.1. 

As r s increases, the PW and Stoll et al. approxima- 
tions tend to the same limit, which is rather different 
from our result. Fig. || suggests that, even if we take a 
conservative approach and fully trust only our r s < 10 
spin-resolved contributions to e c , the common PW and 
Stoll et al. low-density tail hardly matches the QMC 
data. 



IX. CONCLUSIONS AND PERSPECTIVES 



We have proposed a new, analytic, spin-resolved, 
static structure factor and pair-correlation function for 
the unpolarized jellium which works in the density range 
r s < 10. Our model functions fulfill a wealth of known 
analytic properties of their exact counterparts, nicely in- 
terpolate the most recent -arid complete QMC data of 
Ortiz, Harris and Ballone,E3 and consistently yield ac- 
curate correlation energies. They can be of interest 
to build up beyond-rLDA. exchange-correlation energy 
density functionals,E3o~E3 for the magnetierfesponse of 
the unpolarized homogeneous electron gas, LEU and also, 



58, for the e 



within the theory developed in Refs. 
correlation in real materials. As a byproduct, we have 
obtained two correlation energy models which work well 
in the entire r s < 00 density range. 

In further developments we plan to extend our pro- 
cedure to the partially polarized jellium and to lower 
densities (r s > 10). A small Fortran code aimed at 
the numerical evaluation of our functions [Eqs. (|27|), 
([ll]), (|Al]), (|A|), (jBl|)] can be obtained upon request 
to Giovanni . BaclieletOromal . inf n . it. 



B. Spin— resolved 
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APPENDIX A: PAIR-CORRELATION FUNCTIONS IN REAL SPACE 



The expressions of Eqs. (|27|) and (|41|) correspond in real space to 



480 



a? 


'10395 




64 


945 


945 



>645 T , 585. r , , 2 705. 45, t , , 4 (a^p) 5 



r [ , 



315, 



64 



(aTi) 3 
+3^a^(-l)" +1 



64 r 16 



345 
~64~ ( 



33 
64 



64 



n=l 
,TT, 



9(6T1) 



n+l 



(6T1)2 



„TT 
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1722(a TT p) 5 - 68(a T V) 6 + (a TT p) 7 



135135 + 135135a TT p - 270270(a TT p) 2 + H4765(a TT p) 3 - 20370(a TT i o) /1 



TT 
(all) 5 



31185 - 31185a TT p + 6930(a IT p) 2 + 17325(a T1 p) 



(Al) 



,TT^3 



-6930(a TT p) 4 + 938(a TT p) 5 - 52(a T V) 6 + (a T V)' 



.11 



14175 + 14175a TT p + 1890(a TT p) 2 + 



-2835(a TT p) 3 - 882(a TT p) 4 + 378(a TT p) 5 - 36(a TT p) 6 + (a U p) 7 



+3]T4T(-i r +i 



Qn+l 

9(6TT)«+i 



(6TT)2 



(A2) 



APPENDIX B: OPTIMAL FIT TO THE QMC CORRELATION ENERGY 

Since the Perdew-WangS (PW2) form is simple and physically motivated, we slightly modify it by introducing ou& 
more free parameter which grants us enough flexibility to accurat ely fi t the new data by Ortiz, Harris and Ballone.113 



We also include the exact r s lnr s and r s coefficients (see Subsec. [ID) in the high-density expansion of the functional 
form, that now reads: 

e c (r s ) = -2A(l + a 1 r s + a 2 r 2 ) In ( 1 + 1 ) . (Bl) 

This modified PW2 form provides a much more drastic separation between the high- and low-density regime with 
respect to the_jariginal PW2 one. Such a separation is crucial to obtain a good fit which both reproduces the new 
QMC energiesEj at the highest densities and avoids undesired effects on the low-density regime (such as a negative 
coefficient for the r s 3 ^ 2 term). The parameters A, /3±, /?2, @3 and ai are fixed by imposing the high-density expansion 
of Eq. @: aip C/A, fa = 0.5/A exp(0.5S/A), /3 2 = 2A(3 2 and /3 3 = 0.5/?i(8/? 2 A 4 -CB + DA)/ A 2 . A best fit to 
new QMC datall3 gives for the 4 free parameters: ai = 5, /?4 = 45, = 32, (3q = 12.7. The resulting low-density 
expansion is -0.39/r s + 0.99/r 3/2 . The maximum absolute error is 1.6 mRy, while the maximum relative error is 
2.4%. 
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FIG. 1. Spin-resolved pair-correlation function of the unpolarized homogenous electron gas plotted against the electron 
separation r scaled by the density parameter r s for eight different values of r s between r s — and r s = 10. Solid line: this 
work; dots: QMC data; dashed line: Perdew-Wang model. 
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FIG. 2. Antiparallel-spin correlation function at zero intereleetronic distance plotted against the density parameter r 3 . Solid 
line: this work; dashed line: electron-electron ladder evaluation.^] In the inset the exact high-density expansioro is also shown 
as dots. 
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FIG. 3. Static structure factor and magnetic structure factor of the unpolarized jellium plotted against the wavevector q 
scaled by the Fermi wavevector qp for the same eight r s values as Fig. |l|. Solid line: this work; dots: QMC data; dashed line: 
Perdew-Wang model. 
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FIG. 4. Total (solid) and spin-resolved (dashed) correlation energy of the unpolarized homogenous electron gas plotted 
against the density parameter r s . QMC data are also shown as dots. 
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FIG. 5. Our parallel-spin contribution to the total correlation energy compared to the Perdew-WangEl and Stoll et al. 
approximations. e c (r s ,£ = 1) is from Ref. uU 
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